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ABSTRACT 



Increasing the severity of the stress history of a structure reduces it's service life. 
Feasibility studies to increase the zero fuel weight of the P-3 Orion depend heavily on the 
resulting decrease in service life of the wing box and airframe. One option of extending the 
service life of existing aircraft is through the replacement or augmentation of critical 
structural members with composite materials. Since structural composites do not yet have 
adequate service life statistics, life predication must be through probability modeling. Such 
modeling can begin with experimental data on accelerated testing of fiber life under several 
sustained load levels. This data can be the basis for an appropriate strength-life model of 
the fiber, which can in turn be related to the strength-life model of the composite by the 
local-load sharing model. The local load sharing model captures the physical failure 
sequence of fiber failure within a composite. Such a strength-life model, when combined 
with structural analysis, can be used to predict an airframe's service life under the changed 
conditions associated with the zero fuel weight increase. 
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I. INTRODUCTION 



The use of fiber reinforced composites for primary aircraft structures has seen 
enormous increase over the past 10 years. In the absence of sufficient field data due to the 
newness, reliability estimates of these components have proven to be difficult to obtain. 
The reliability of a structure is the probability that it will survive a given load history. 
Without a firm reliability model the predicted service life of a composite structure will not 
be accurate. This problem may result in both low aircraft availability and a degradation of 
aircraft flight envelope to maintain aircraft safety. If an increase in airframe loading is 
desired, then the effects on service life must be estimated beforehand. What is needed is a 
method for accurately predicting the reliability of a composite structure for a variety of load 
histories. 

Traditional methods of predicting the service life of mostly aluminum aircraft 
structures have been based on over sixty years of historical laboratory and field data. 
Statistics of how a material survives a life cycle of loading allows for construction of S-N 
type curves from which service life is predicted. Such experienced based relations between 
strength and an equivalent life cycle is thus known for many metals. For advanced 
composites, the historical data base reaches back approximately twenty years. However, 
advanced composite materials with new fiber and matrices are continually being developed. 
The field data do not belong to the same set and cannot be merged to formulate empirical 
life prediction models The relationship between the strength and life characteristics of these 
composites remains unknown. 

One method of filling this void in composite statistics is through distilling the salient 
failure processes into a mathematically tractable probability model. This approach is 
general only when the failure mechanism for the basic material element remains unchanged 
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for different composites. The largest such basic element appears at the fiber level. The 
failure of a single filament fiber is a homogenous Poisson, process which can be 
adequately modeled by a weakest link Weibull distribution. However, the failure process 
at the composite level is not homogeneous. Under tensile conditions, it is sequential in 
stress level or in elapsed time starting with failure of weak fibers in spatially dispersed 
sites. The local failures are isolated by the load transfer action of the matrix binder. As 
stress or time increases, the number of failure sites increases until chanced clustering 
ultimately leads to global catastrophic failure. This failure process has been represented by 
the local load sharing model. The relationship in strength between a composite and the 
parent fiber has been verified experimentally. It has been demonstrated that the strength 
statistics of the fiber can predict the probability of failure of the composite under increasing 
stress. This investigation explores a parallel approach to the modeling of composite failure 
under increasing time. The statistics of fiber failure in time are being characterized by 
stress rupture experiments. The resulting fiber life statistics will then be incorporated into a 
load sharing model with time as a random variable to predict the life of the composite. This 
strength-life model for the composite can then be used to project the life of the structure. 

Manufacturing considerations may require the use of new material before little or no 
data testing is performed. Because of the lead time required to obtain strength-life data, an 
accelerated strength-life experiment must be performed. This involves the use of higher 
than expected load levels to produce more timely failures. If the failure mechanism for the 
test specimen is homogeneous, then the accelerated testing will produce results similar to 
those experienced during service life. 

The difficulty in performing such a task is apparent when one realizes that a certain 
dichotomy exists between testing of strength and life. That is, once a fiber is tested in 
strength it is physically impossible to retest that fiber in life. This seemingly trivial 
observation makes the direct statistical correlation between the strength and life of fibers 
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from the same population a true challenge. This correlation is the intent of an ongoing 
research experiment at the Naval Postgraduate School (NPS) Advanced Composites 
Laboratory. Two graphite fiber populations are simultaneously being tested to obtain 
strength and life statistics. Potential strength-life models will then be tested against these 
statistics. Once a strength-life correlation for fibers has been determined, the same relation 
will be investigated using composite strands from the same population. This research will 
verify the strength-life degradation relation between a composite and the parent fibers from 
which it was made. 

The intent of this investigation is to present a methodology to correlate fiber strength 
and life statistics initially with a currently available model and explore an application to 
service life prediction. Laboratory data from both strength and life experiments will be 
used to formulate mathematically consistent degradation models. These models must then 
be extended to account for a broad range of load histories. Once the relationship between 
life and strength is understood, a probabilistic estimate of the life of a fiber for a specified 
load history may be determined. To illustrate the usefulness of such a methodology, a 
current Navy operational problem will be analyzed. A proposed increase in the zero-fuel- 
weight of the P-3 Orion requires addressing the affect this will have on service life. 
Replacement or augmentation of critical structural areas of the airframe with composite 
materials is one possible solution. The ability to assess what affect this will have on the 
airframe service life is indeed most desirable. An analysis of this problem is included as 
Appendix B. The main body of this investigation looks into the theory and experimental 
requirements for predicting the service life of composite materials. 
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II. BACKGROUND 



Several topics will be presented as background information for analyzing strength and 
life data and constructing appropriate models. The failure mechanisms of a fiber and a 
composite must first be understood from both the strength and life viewpoint. This will 
allow for selection of an appropriate model which describes the failure, or breakdown 
process over time. Life testing of composite fibers requires a knowledge of the affects of 
data censoring in an accelerated life test. These concepts will assist in the formulation of 
appropriate experiments and an accurate strength-life dependency model. First, the nature 
of the P-3 zero fuel weight problem will be analyzed further. 

A. WHY SERVICE LIFE PREDICTION? 

As the fleet of aircraft known as Naval Aviation continues to age with little prospect 
for replacement in the near future, many aircraft structural components are approaching the 
end of their original designed service lives. Replacement with replica components is 
always a possibility, however, years of adding weight in the form of new armament and 
avionics systems have made original design considerations obsolete. Such is the case of 
the P-3 Orion, where calls for increasing the zero fuel weight by close to 18,000 pounds 
raise serious questions about the effect this will have on structural integrity and service 
life. 

The zero fuel weight (ZFW) of an aircraft is comprised of the aircraft structure and 
payload. This weight is used to determine maximum weight for landing, when fuel 
reserves have been depleted to a safe level. Increasing the ZFW of an aircraft generally has 
the effect of increasing the landing weight if fuel reserves remain constant. In the case of 
the P-3, a ZFW subjects the wing box in particular to a large increase in load. This specific 
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problem was examined in by Lt. Culpepper [1], who analyzed the effects of overloading 
which accompanied increasing the P-3 ZFW. 

A consequence of this overload is an expected life shortening of internal members of 
the wing structure. Solutions to this problem range from augmentation of current 
structures with like materials to replacement of the wing box with composites. The area 
under consideration is shown in Figure 2.1. Regardless of the method of strengthening the 
wing box, their effect on the aircraft service life must be estimated by known methods of 
probability analysis. Figure 2.2A illustrates conceptually that the service life of the P-3 
fleet at current ZFW may be described by a cumulative distribution function (CDF). This 
curve indicates that at some time tj, the probability that a certain percentage of the fleet is 
grounded is 10%. If the ZFW is increased, a higher stress will produce a shift of the CDF 
to the left. The time for grounding 10% of the fleet is also shifted left to some time tj' as 
shown in Figure 2.2B. The difference between tj' and ti is the reduction in service life at 
the 10% retirement level. If the dependency of the CDF on stress level is known, then fleet 
availability at any service stress level (or the ZFW) can be projected. 

This is illustrated by a three dimensional representation of the CDF dependency on the 
ZFW in Figure 2.3. The projection of the family of CDF curves on a plane gives rise to the 
equivalent of the widely used S-N curve. Determining the CDF curves for a range of stress 
values would produce a surface of the form found in Figure 2.4. Knowing the shape of 
this surface, one could postulate the effect of changes in aircraft loads on service life, a 
most desirable tool indeed. 

The probability of failure surface as depicted by Figures 2.3 and 2.4, can be produced 
by statistical data for material failure under various load histories. For aluminum structures 
this is not a problem, as the failure mechanisms for aluminum are well known and are used 
to produce the S-N curves found in design handbooks. Composites do not enjoy such 
status, for the availability of composite failure statistics is based on limited service history. 
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The surface in Figure 2.4 is completely hypothetical for a composite material since the 
relation between strength and life has not been statistically verified. 

As applied to the problem at hand, if the P-3 wing box is structurally augmented with 
aluminum, then the effect on the life can be estimated. However, if composites are used, a 
means of determining a CDF surface must be devised. The next two sections cover more 
required background material before the experimental and modeling means for producing 
the CDF surface are discussed. 

B . COMPOSITE FAILURE MECHANISMS 

In order to understand the relationship of composite failure in both strength and life it 
is important to understand the process of failure. A fiber reinforced composite material is 
comprised of two principle components, the fiber and matrix. The fiber supports the load 
applied to the composite structure, while the matrix acts in shear to transfer loads once 
carried by broken fibers. This action provides the micro-redundancy associated with the 
reliability of a composite. These concepts are fundamental to describing a composite in 
failure and will be discussed henceforth. 

1 . Composite Strength Model 

In order to describe the failure of a fiber composite one must first look at the 
failure of an individual fiber. A fiber may be described by the "weakest link" model, in 
which an item fails when one of a series of integral components fail. Thus a fiber may be 
thought of as a series of links, where the failure of one link will cause the failure of the 
entire fiber. This well documented failure mechanism is characterized by the Weibull form 



where a and p are shape and location parameters respectively. If these parameters are 
known, the probability of fiber failure for a given stress level may be determined. 




( 2 . 1 ) 
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FIGURE 2.1 P-3 WING BOX 
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FIGURE 2.2 EFFECT OF ZFW INCREASE ON SERVICE LIFE 
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FIGURE 2.3 STRESS-TIME SURFACE CONSTRUCTION 

Failure of a fiber occurs when the stress in any of these hypothesized links is 
greater than the strength characterized by Eq. 2.1. When combined with a matrix the effect 
of the failure is isolated and the load once carried by the broken fiber is transferred to 
surrounding fibers through shear. This concept, termed local load sharing, was initially 
developed by Rosen [2], who quantified the distance a broken fiber was unable to carry a 
load. This distance, termed the ineffective length 6, is illustrated in Figure 2.5. When a 
fiber breaks, the axial stress o in the fiber abruptly goes to zero while the shear stress x in 
the matrix maximizes at a similar rate and serves to transfer the load to adjacent fibers. 
Initially, fiber failure sites are dispersed throughout the composite as reflected by the fiber 
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strength statistics. As the applied stress increases, weaker fibers continue to break with 
their load transferred to neighboring fibers. The number of fiber breaks will eventually 
produce a clustering of fiber failure sites. When the ineffective lengths begin to overlap the 
load sharing mechanism will cause further breaks until catastrophic failure of the composite 
occurs. Thus, the local load sharing model describes the micro-redundancy of a composite 
in failure. 
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FIGURE 1.5 TENSILE FAILURE MODEL [2] 

The relationship between the strength of a fiber and composite can be 
characterized by the fiber strength statistics and the matrix effectiveness. Harlow and 
Phoenix [3,4] have described this relationship in terms of a probabilistic model based on 
the number of failed elements k adjacent to a fiber. This model is known as the Harlow 
and Phoenix Local Load Sharing Model. As k increases, the strength of the adjacent fiber 
will eventually be reached, resulting in failure. This process will lead to the catastrophic 
failure sequence described above. Harlow and Phoenix perform a recursion analysis for 
increasing values of k and provide numerical results which may be used as a benchmark 
for further study. 

Experimental verification of the fiber-composite relation was performed by 
Johnson [5] using strength data from Hercules Magnamite AS-4 graphite fiber and 
composite strands. Johnson concludes that the Local Load Sharing Model provides a 
means for predicting composite structural reliability using fiber strength statistics. A visual 
representation of fiber-composite relation may be seen in Figure 2.6. The data denoted by 
the larger circles represent fiber strength data and the smaller circles represent composite 
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strength data. The solid line is based on the Local Load Sharing Model, which clearly 
describes the failure characteristics of the experimental composite data. 

2 . Composite Life Model 

The life prediction of a composite is based on the time dependency of the 
parameters in Eq. 2.1. This idea was initially discussed by B. Coleman [6] in his theory 
on the time dependence of the mechanical breakdown of fibers. Coleman developed the 
notion of a breakdown function tc(o(t)) to describe the degradation of a fiber under a 
specified load history. He concluded that the failure rate, or hazard l P is a function of the 
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breakdown 'P(k). This generalized form was used by Phoenix and Wu [7] to characterize 
the life of a Kevlar-49/epoxy composite and is discussed in detail below. 

When the failure mechanism is homogeneous, the hazard ^(t) is characterized by 
a power form of the Weibull distribution, W(x)“T a . Where t is the reduced time as 
dictated by the break down rule K and the stress history o(t). This reduced time is a 
fractional measure of a fiber's life as degraded by the severity of the stress history and is 
cumulated in the integral form 



x - 



|/\{o(')}dt 



(2.3) 



t is a normalizing constant interpreted as the intrinsic life with the dimensions of time. 

Once the reduced time has been convoluted into the breakdown rule, the joint 
distribution between stress and time is described by: 



F(t|o) - 1 - expj-W ji K[o(t)]dt|j 



(2.4) 



For this formulation, the breakdown rule the only unknown. The breakdown rule may be 
of any form which can be reconciled with observed data. Two simplest possible forms are 
the power law and exponential law. The power law is based on a multiplicative, or flaw 
growth, failure mechanism of the form 



K 




(2.5) 



Whereas the exponential law is based on activation kinetics, or a flow growth based failure 
mechanism of the form 




( 2 . 6 ) 



Definite validation of the proper breakdown rule requires a data base consisting of wide 
ranges of stress history. For the purposes of this investigation, the power law, Eq. 2.5, 
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will be used without loss of generality. The intent is to present a methodology for 
validating the appropriateness of a chosen form of the breakdown rule. 

Substitution of Eq. 2.5 into Eq. 2.4, an integration for a constant stress history, 
produces a relation of the form 



Note that a variety of load histories o(t) can be represented by piecewise combinations of 
Eq. 2.7 and 2.8 for specific time intervals. Complete derivations of Eq. 2.7 and 2.8, along 
with the exponential form , Eq. 2.6, are provided in Appendix A. 

C . ACCELERATED LIFE TESTING 

The theoretical relationship between a composite and a fiber in strength has now been 
established. Two well known forms of a strength-life relationship have also been derived 
and discussed. What is needed now is a manner in which to utilize these concepts to 
correlate experimental results, including both strength data and life data. In order to test 
any model, a sample from a single statistical population must be tested in both strength and 
life. The amount of time and tedious management required to obtain a statistically adequate 
set of life data is enormous. A method of decreasing the time, while maintaining acceptable 
confidence limits, is needed. 

1 . Data Censoring in a Life Test 

When performing a life experiment, very seldom is it the case that all experimental 
samples are realized in failure time. Data from life tests frequently have from three to four 
magnitudes of scatter, so time limitations dictate that the termination (or censoring) of long 
life samples be completed. Censoring is used to terminate an experiment once enough 




(2.7) 



For a simple load history of constant load rate, L, the relation of the parameters are 




( 2 . 8 ) 
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information is obtained. It is also used to obtain partial information from data which is 
tainted so that the experimental effort is not completely lost. In order to obtain data in a 
more timely fashion, over stressing the test specimen above service loads will further 
accelerate the test. As an experiment progresses, updated estimates of the parameters of 
interest may be completed until their confidence is deemed acceptable. If reliability is of 
major concern, the characteristics of initial failures (lower tail distribution) are of more 
interest, negating the usefulness in continuing a life experiment. 

2 . Types of Data Censoring 

Censoring occurs in two basic forms, right censoring and left censoring. If a 
specimen is known only to have survived longer than the censor time, the data is said to be 
censored on the right . Similarly, if a specimen is known to have failed before some 
known time, the data is censored on the left. A data set with a combination of the two is 
said to be interval censored. 

The proper time to terminate a right censored experiment may be thought of in two 
ways. First, an experiment may be terminated after a predetermined time has elapsed. This 
method is known as time censoring, or a Type I censor. Secondly, an experiment may be 
terminated after a certain number of failures have occurred, keeping the time to the fixed 
number of failures random. Such data is known to be failure, or Type II, censored [8]. If 
a specimen under test is failed prematurely through some outside influence, the data is said 
to be a Type I-v censor. 

3 . Parameter Estimation 

The primary goal of conducting a life experiment is to estimate the distribution 
parameters which describe the sample. In the case of the Weibull distribution we want to 
establish the shape parameter a and the location parameter of Eq. 2.1. One method of 
determining the most appropriate parameters is through the use of Maximum Likelihood 
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Estimators. The probability, or likelihood that a selected a and P are in fact the correct 
parameters is described by 



"IWflM'#*' 1 ' 



(2.9) 



Where m is the number of exact, or realized data, r is the number of right censored data, 
and / is the number of left censored data. The sample size is defined as n=m+r+l. An 
expansion of this equation for the Weibull distribution is included in Appendix C. 

The a and (3 which have the highest probability of being the true parameters are 
called the maximum likelihood estimators and may be obtained through the simultaneous 
solution of the following equations: 



dL 

da 

3L 

dp 



0 

0 



( 2 . 10 ) 



Depending on the form of Eq. 2.9, Eq. 2.10 may be solved analytically. If this is not the 
case, numerical methods must be used. One method of circumventing this difficulty is 
through graphical interpretation of the likelihood surface. Given a data set x lt and an 
appropriate range of a and P, Eq. 2.9 may be solved for each combination of a and p. 
Plotting a 3-D surface, or contour, the maximum may be readily seen in the form of a peak. 
This effect is demonstrated in Figure 2.7, where the likelihood for a set of N=64 simulated 
points with a=5 and p=20 is plotted. These parameters are typical of the fiber strength 
data displayed in Figure 2.6. 
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Likelihood Contour Levels 




beta (normalized) 



Likelihood Surface 




FIGURE 2.7 LIKELIHOOD CONTOUR AND SURFACE 
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4. Effect of Censoring on Parameter Confidence 

A primary consideration in determining the MLE is the effect that censoring will 
have on accuracy. J. Coleman [9] investigated the effect of censoring on an experiment 
and the optimization of information from varying censor time. He concluded that a 
recursive censoring scheme would optimize the information gained from experiments with 
high variability, such as a life experiment. Hence, more information may be gained from 
censoring multiple experiments instead of allowing a single experiment to continue until 
completion. 

With that in mind, Figure 2.8 demonstrates the effect of censoring the same 
simulated data used for Figure 2.7. It can be seen that as the number of points censored on 
the right is increased, the estimated parameters lose their accuracy. However, the peak of 
the likelihood surface is still seen to contain the underlying parameters. If the same 
censoring scheme is performed on data with higher variability, a decrease in parameter 
confidence occurs. In Figure 2.9 the same underlying rank for the data set in Figure 2.8 is 
used to produce a data set where a«.2, more typical of a life experiment. The result is a 
noticeable effect of right censoring data on the estimated parameters. If in addition to the 
right censor a left censor is employed, the confidence of the parameters also decreases, but 
not as dramatically as the right censor case. This effect is displayed in Figure 2.10. The 
low confidence levels associated with shape parameters less than one make the 
determination of the MLE quite difficult. Simulations using expected parameter values 
allow the evaluator to recognize subtleties in the likelihood surface caused by censoring. 
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A. Likelihood Contour and Surface: 30 Points Realized 





beta (normalized) 

B. Likelihood Contour and Surface: 20 Points Realized 





beta (normalized) 

C. Likelihood Contour and Surface: 10 Points Realized 





1 1.5 

beta (normalized) 



FIGURE 2.8 EFFECT OF RIGHT CENSOR, N=64, a=5 
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A. Likelihood Contour and Surface: 30 Points Realized 





beta (normalized) 

B. Likelihood Contour and Surface: 20 Points Realized 





beta (normalized) 

C. Likelihood Contour and Surface: 10 Points Realized 




beta (normalized) 




FIGURE 2.9 EFFECT OF RIGHT CENSOR, N=64, a-.2 
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A. Likelihood Contour and Surface: 5 Points Left Censor, 30 Points Right Censor 






beta (normalized) 

C. Likelihood Contour and Surface: 20 Points Left Censor, 30 Points Right Censor 





beta 



FIGURE 2.10 EFFECT OF LEFT AND RIGHT CENSOR, N=64, a=.2 
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III. FIBER STRENGTH-LIFE EXPERIMENT 



Ongoing research at the NPS Advanced Composites Laboratory has produced both 
strength and life data for two fiber populations made from the same manufacturing process. 
The intent is to produce a statistical data set which will validate or refute any proposed 
strength-life degradation model for both fiber and composite. Strength tests have been 
completed on both fibers and composite strands with results documented in References 8 
and 9. Continued efforts will produce additional data from life tests in progress. 

A. TEST SPECIMENS 

Fiber and composite specimens were produced from two production spools of AS-4 
Hercules Magnamite high strength graphite, designated 008 and 019. These two spools 
were treated as two distinct populations for both the strength and life tests. Composite 
strands were constructed from bundles of approximately 3000 fibers and strength tested in 
gauge lengths of 2 and 10 inches. Composite strand life tests have yet to be completed. 
Fiber samples were strength tested individually using gauge lengths of 50 mm. Fiber 
samples being tested in life use the same gauge length. 

B. STRENGTH TEST RESULTS 

Reference [10] contains documentation for strength tests of fibers from the 008 and 
019 AS-4 spools. Over 300 specimens were tested in strength and analyzed to determine 
the MLE for both spools. Parameter estimations are summarized in Table 3.1. Reference 
[5] documents the strength testing of 82 composite strands. Similar MLE calculations were 
performed and are also listed in Table 3.1. The subscript "r" denotes that the parameters 
are for strength tests conducted under a constant loading rate, or ramp load condition. 
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Comparisons were made between the fiber and composite strength data, verifying that if 
fiber strength statistics are known, a composite probabilistic strength is obtainable. 

TABLE 3.1: STRENGTH TEST MLE FOR AS-4 SPOOLS 008 AND 019 [8,9] 



AS-4 Spool/test 


A 




008 fiber 


5.05 


18.4 


019 fiber 


5.17 


17.4 


008 10" strand 


15.83 


54.53 


008 2" strand 


33.62 


57.07 


019 10" strand 


17.48 


48.34 


019 2" strand 


14.75 


56.52 



C . LIFE TEST 

Fiber life testing began in August, 1991 and continues to this day. Original design 
considerations placed eight groups of 64 fibers at four different constant load levels. This 
arrangement resulted in a total of 512 fibers on test with half from each spool. This 
arrangement will produce two statistically independent data sets which are relatable to 
completed strength tests. Individual AS-4 fibers were suspended from low compliance 
springs with a tensile load produced from weighted vials. Four different load levels were 
chosen with consideration given to future calculations of the degradation relation. These 
levels were taken to be 14.17 gm., 13.15 gm., 12.20 gm., and 10.96 gm. Fibers were 
loaded using a calibrated load cell as a platform to measure failures on loading. During the 
initial loading cycle a number of fibers were tested in strength to verify the statistical link 
between the strength and life data [11]. 

Fiber failure times are monitored by a data acquisition system controlled by a PC-AT 
computer. Failures are sensed by the closure of an infrared switch surrounding the 
compliance spring. Fiber life duration's are then computed and compared with written 
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records prior to addition to the fiber life data base. Reference 11 describes the original 
loading scheme in more detail. A complete description of the loading apparatus with fiber 
loading procedures is included in Appendix D. Instructions for use of the monitoring 
system may be found in Appendix E. 

Results from the fiber life experiment are still in their infancy. With the exception of 
the 14.17 gm. load level, current realizations in life are too few to produce accurate 
parameters. Initial fiber loading included two additional rows of 14.17 gm loaded fibers. 
The intent being to produce additional data at the higher load to improve confidence of 
initial parameter estimation. This should provide a better estimation of the parameter p, to 
be discussed in Chapter IV. Subsequent loadings at the lower stress levels have yet to 
produce substantial realizations. 
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IV. STRENGTH-LIFE CORRELATION OF COMPOSITE 

FIBERS 



Probability Modeling of a composite structure is a two-tiered process. First, statistical 
data of both strength and life characteristics of the composite must be gathered. Second, a 
mathematical model of the strength life relationship must be developed. As discussed in 
Chapter III, research conducted on composite and fiber from the same statistical population 
has validated a relation in strength between the two. A correspondence between strength 
and life characteristics is being investigated. Since the failure mechanism for a fiber is 
much simpler than a composite, the use of a more explicit fiber model will be investigated. 

A. STRENGTH-LIFE TEST CONCEPTS 

A fiber is said to be tested in strength if stress at failure is the random variable, while a 
fiber is tested in life if the time to fail is the random variable. Before developing a strength- 
life model, the relation between the two tests must be understood. Both the strength and 
life tests may be viewed in either a deterministic or probabilistic manner. The stochastic 
nature of fiber failure dictates that the probabilistic view is appropriate. For a sample of 
size N, a distribution will evolve between the weakest and strongest specimens. Figure 4.1 
depicts this for both the strength and life tests. 

A simple load history for a strength test, for example, is the application of a constant 
load rate, or ramp load until failure. Failure load is measured for a sample of size N and 
statistical analysis provides the appropriate descriptive parameters. In this case, the failure 
stress depends on the loading rate. Likewise, a simple life test may be performed by 
loading a specimen to a predetermined stress level, which is held constant until failure. 
Depending on the target load level, there may be a number of specimens which do not 
survive the loading procedure. These specimens are said to be realized in strength and left 
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censored with respect to the life test. The number of specimens realized in strength will 
depend on both the loading rate and the stress level chosen for the life test. Example load 
histories are depicted in Figure 4.1. 

To investigate the underlying relationship between strength and life, the tested samples 
must belong to the same populations. In an idealized thought experiment, the same sample 
could be tested once in strength and then again in life. If this were physically possible, the 
relation between the two tests would be apparent with the slope of degradation obtainable 
from one single fiber sample. Figure 4.2 illustrates this concept for a sample of six 
idealized fiber specimens. Note that under this thought experiment the relation between the 
fibers in the sample remains unchanged through the various load histories. Since one 
sample cannot be tested more than once, the deterministic strength-life degradation cannot 
be observed. To better visualize this concept, the simulated sample from Chapter II with 
N=64, a«5, and {W20 was tested by first incrementing stress rupture load levels and then 
incrementing loading rates. The resulting CDF surfaces were plotted with the results 
shown in Figure 4.3. It is readily apparent that the different load histories produce 
different surfaces and hence different probabilities of failure. In other words, the so called 
'S-N' curve is load history dependent. 

To simplify analysis and allow for evaluation of the strength-life degradation, all load 
histories will be reduced to the fundamental building block of a stress rupture load history. 
A stress rupture experiment is a test in which a sample is loaded at an infinite loading rate 
up to a target load level and sustained until failure occurs. Since in an actual stress-rupture 
experiment a finite loading rate is used, the samples have a combined history. In fact, 
some samples fail during a constant rate loading phase, while some fail during the constant 
load phase. To simplify analysis and allow for evaluation of the strength-life degradation, 
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all load histories will be reduced to the fundamental building block of a stress rupture load 
history. A stress rupture experiment is a test in which a sample is loaded at an infinite 
loading rate up to a target load level and sustained until failure occurs. Since the strength 
and life tests of Figure 4.1 contain a constant loading rate as the entire or part of their 
history, they are not considered failures in stress rupture. Since the strength and life tests 
of Figure 4.1 contain a constant loading rate as the entire or part of their history, they are 
not considered failures in stress rupture. 

Obviously, loading a fiber instantaneously is not physically possible. To allow for 
this, a mathematical transformation must be performed to set At=0 as represented in Figure 
4.4. For the strength test case an arbitrary time t may be chosen to be the intrinsic time at 
which the loading occurred. This is an important concept since the logarithmic time scale 
contains no time zero. Here the intrinsic time is interpreted as the stress rupture time at 
which a strength is defined. The larger the designated t , the more time the specimen has 
been under stress, and the weaker the sample has become. This assumes that there is a 
degradation of strength over time, i.e. the sample does not become stronger over time. 
Therefore, an appropriate model must demonstrate that an intrinsic time may be chosen 
arbitrarily. If the chosen form of the breakdown rule is linear, then the intrinsic strength is 
the intercept value at the desired intrinsic time. 

For the life test case however, the transformation of specimens which failed on 
loading to the stress rupture domain requires a knowledge of the strength-life relation. For 
the loading failures to be considered part of the life test, an intrinsic time must be chosen 
such that the specimen failed at the precise moment of reaching the life test load level. This 
transformation will merge the data which failed during loading to the stress rupture data, 
thereby alleviating any censoring on the left. Removing left censoring should produce 
higher confidence in the estimated parameters. 
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When predicting the service life of a composite a variety of load histories must be 
considered. For analysis purposes it is desired to relate these operational load histories to 
the basic stress rupture form. For example, a typical flight load history for the wing box of 
a P-3 Orion may be depicted as shown in Figure B.l. What is desired is an equivalent life 
for a constant stress level (stress rupture), noting that an increased stress level shortens the 
equivalent time. This type of load history may be broken up into segments of constant 
stress and constant load rate which may be individually transformed to an equivalent stress 
rupture history. A detailed example of this process is provided in Appendix B. 
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Stress Rupture Surface 




Constant Load Rate Surface 




FIGURE 4.3 EFFECT OF LOAD RATE ON STRESS-TIME SURFACE 



30 




FIGURE 4.4 DATA MAPPING TO STRESS RUPTURE DOMAIN 
B . DEVELOPING THE STRENGTH LIFE MODEL 



Given the appropriate break down rule k and it's parameters, the expected life for any 
load history can be estimated. The fiber data realized in strength and those realized in life 
can be used to establish the function form of the break down rule K and provide an 
estimation of the parameters in k. Since actual life data will not be completely realized in 
the near future, simulation is used to establish analysis procedures and assess the 
usefulness of various forms K. As more data become realized, insight is gained into the 
strength-life relationship, allowing for modifications of the breakdown rule. The Monte 
Carlo simulation is used to simulate the completed life experiment. Censoring and analysis 
techniques are then evaluated so that interpretation procedures may be improved while the 
experiment is still ongoing. 
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Two assumptions have been made to focus analysis of the data and document the 
procedure. First, the failure potential developed by B. Coleman [6] is assumed to be of the 
correct form. Second, our analysis of the breakdown rule will center around one of two 
forms developed, either the power law or the exponential law. The form chosen should be 
based on the known micro mechanical behavior of graphite fibers during stress rupture. In 
general, failure mechanisms may be placed into two categories: flow growth and flaw 
growth. Flow growth is based on kinetic fracture theory and may be described by an 
exponential form. Flaw growth on the other hand is probabilistic in nature and is modeled 
using a joint probability form such as the power law. Graphite fiber is known to be more 
brittle in fracture, so the power law should be the appropriate model form of the breakdown 
rule. When performing actual data analysis, other models may be investigated. The power 
law will be specifically investigated here. 

After model selection, correct estimations of the model parameters are needed. For the 
case of the power form of the breakdown rule, the unknown parameters are p and A. 
Placing Equation 2.7 in the linear logarithmic form log(y)=mlog(x)+b produces Equation 
4.1. 



log(o) - - -log( + log(A) 

P ' * > (4.1) 

The parameter p is seen to represent the negative inverse of the slope, while A represents 

the ordinate intercept. Note that this interpretation of the parameters was conducted using 

the stress rupture load history. If a different load history is used, the change in domain 

requires a change in interpretation. In this case, a transformation of load history to the 

stress rupture domain will allow for a relatively straight forward analysis of the parameters. 
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C . MODEL VALIDATION THROUGH SIMULATION 

To assist in validation of a proposed break down rule, a software package has been 
developed and is included as Appendix F. Program STRENGTHLIFE allows the user to 
input the form of the breakdown rule as a subroutine and obtain a graphical display in the 
stress rupture domain. A three dimensional output of the cumulative failures surface in the 
stress-time space is also produced for the proposed breakdown rule. If the analysis is 
performed correctly, the simulated value of p should be recovered as data is transformed to 
the stress rupture domain. 

To validate the software, a simulated data set was produced with parameters shown in 
Table 4.1. Shape and location parameters used were based on results described in Table 
3.1. The value of p chosen was based on results found in [7] for Kevlar-49. The 
simulation produced a ramp loading strength data set which was transformed to a four 
stress level life test. These input values closely follow the fiber life experiment in progress 
at the NPS Advanced Composites Laboratory. A fifth, much lower load level, was chosen 
to represent a possible service load history. If the analysis is performed correctly, a linear 
relation with slope p between data points will exist, as shown in Figure 4.2. 

TABLE 4.1 INPUT TEST PARAMETERS FOR DATA SIMULATION 
a 5 

P 20 

p 20 

L -8 gm/sec 

t 10 -15 , 10°, 10 24 

Mapping load levels 14, 13, 12, 11, 2 gm 

As previously mentioned, the intrinsic time may be arbitrarily chosen without loss of 
generality. For the current problem three different values of t were chosen. The results, 
as shown in Figure 4.5, demonstrate that the values for the five sets of life data remain 
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constant as t is shifted. Consequently, the slope remains constant and is taken to be a 
correct value for p. The values of the intercept A, which may be interpreted as intrinsic 
strength, is seen to decrease with time. This is consistent with the assumption of a 
decrease in strength over time. This also demonstrates that the faster a specimen is tested in 
strength, the higher the failure load. Thus a strength test conducted infinitely fast will be 
met with a failure stress level infinitely large. Not a very intuitive thought at all! 

The strength-life relationship demonstrates that if the power law is a correct form, 
experimental strength data may be used to predict failure times at various stress levels. The 
representative service life loading may be adjusted as mission and airframe requirements 
dictate. Since the value for p was simulated, the data transformation represents a simple 
linear mapping. An experimental value for p is required to test the validity of the preceding 
analysis. 




FIGURE 4.5A SIMULATED MAPPING FOR 1=10 ^ sec 
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FIGURE 4.5B SIMULATED MAPPING FOR t=l sec 




log(t) - sec 



FIGURE 4.5C SIMULATED MAPPING FOR t =1024 sec 
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D. ANALYZING EXPERIMENTAL DATA 

The experiment outlined in Chapter III is designed to identify the form of the break 
down rule model. Since the strength and life tests fibers are from the same statistical data 
set, the distributions should be related. Englebert [10] calculated the parameters for the 
fiber strength test. Currently, not enough realizations exists from the life experiment to 
conduct the same analysis. The different types of censoring present in the life experiment 
reduce the accuracy of parameter estimation. A method for minimizing the number of 
censored data would produce the most likely estimated parameters. 

The analysis objective is to obtain the best estimation for the distribution parameters 
using the most information obtainable. To accomplish this, an iterative method, which 
maps all data to the stress rupture domain, would produce the fewest number of censor 
points. Starting with the most basic case, as shown in Figure 4.6A, the data set is left 
censored for fibers which failed on loading and right censored in time for fibers which have 
not failed. This case produces the fewest number of exact data points for calculation of 
MLE. Figure 4.6B demonstrates that if the fibers which failed on loading were transferred 
to the stress rupture domain, they would in effect become exact data points which would 
increase the accuracy of the MLE. If the fiber were to undergo a strength test some time 
during the life experiment, the number of right censored data may be reduced. Figure 4.6C 
indicates that failures outside the stress rupture domain could be transferred back to 
produce additional exact data points on the right. This method, known as proof testing, is 
a technique used to further accelerate a life test and optimize information gained. 

Mapping data to the stress rupture domain is an excellent means for increasing the 
number of exact data points for calculating MLE. However, as mentioned above the value 
of p must be known. To remedy this problem, an iterative approach is used to estimate p 
for the purpose of mapping. Finding the MLE for the censoring scheme in Figure 4.6A 
permits construction of the plot shown in Figure 4.7. A value of p may now be estimated 
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and used to conduct the data transfer of Figures 4.6B and 4.6C. Recalculation of the MLE 
using the greater number of exact data will produce a refined values for a t) (5 t , and p. This 
value of p may then be used to evaluate the appropriateness of the model by mapping (5 S to 
the experimental stress levels as depicted in Figure 4.8. If the values of experimental p t 
and transferred |5 t agree, then the model is validated. If the model is not validated, then 
modifications may be made and the process repeated. Once a model has been validated, 



fiber life predictions based on fiber strength statistics may be produced. 
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FIGURE 4.7 FORMULATION OF p USING EXPERIMENTAL DATA 
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V. CONCLUSIONS AND RECOMMENDATIONS 



The accurate service life prediction of composite aircraft structures is a technology void 
which must soon be filled. Unfortunately all the money in the world cannot buy the time 
required to perform meaningful life experiments. The research being conducted at the NPS 
Advanced Composites Laboratory will provide a data base to formulate a strength-life 
model. Such a model would be the basis for performing accelerated stress experiments to 
reduce testing time. Once an appropriate strength-life model is verified, composite service 
life may be estimated using fiber strength statistics. 

The following topics are recommended for follow on research: 

1 . Determination of accurate MLE for the AS-4 life test fibers such that a strength-life 
slope of degradation may be calculated. 

2. Study the effects of various forms of the break down rule to identify the appropriate 
strength-life model. 

3 . Initiate the testing of AS-4 composite strands in life such that a composite-fiber 
relation in life may be verified. 
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APPENDIX A: INTERPRETATION OF THE BREAKDOWN 

RULE 



The breakdown rule k(o), as discussed in Chapter II, must capture the principle failure 
process. The challenge is to develop a tractable form which adequately depicts the 
degradation of life under a range of load histories. The power law and the exponential law 
are two such candidates. Repeating for clarity, the power law is of the form 

p 

(A.1) 

and the exponential law, 







(A.2) 



As a starting point we apply this failure mechanism to the simple load history of stress 
rupture. The sample is loaded at an arbitrary time to to a constant tensile load oj. The 
resultant degradation in life may then be investigated using the power law and exponential 
law breakdown rules as discussed. Figure A-l depicts this loading scheme. 




FIGURE A-l STRESS RUPTURE LOADING 
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Using the power law breakdown rule 



T - 7jM°w} d ' - iC(f) l>d ' - {jfr where '• ■ 0 

Substituting this into the hazard function we obtain a reliability of, 

»P / ( ' a 



R(/) - exp[-W(/)] - exp 



-(f) (*) 



By comparison with an exponential distribution of the form 



exp 



-fe) 



. we find that a=at and 






(A-3) 



(A-4) 



(A.5) 



which may be represented in a logarithmic manner for ease of plotting as, 



log(o) - --log(&-\ + log(A) (A.6) 

p V t } 

a linear relation in a log-stress, log-time plot. Values for the location parameter pt and 
intrinsic life t are obtained from a life experiment. For a given stress load we are left with 
p and A as unknown parameters which, when combined graphically with different loads, 
may be obtained from the slope and intercept information. 

Similarly, using the Exponential Law Breakdown Rule one obtains, 



x - jjM o w} d ' - rX'Wf ) d ' ■ i exp (§) t where '» ■ 0 (A - 7) 

which produces a reliability of 



R(t) « exp 




comparisons with the exponential distribution again show a=at and 



(A.8) 
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(A.9) 



This may be represented in a logarithmic manner as, 



o 



--Clogj'y + log(B) 



(A.10) 



which is linear in a stress vs. log-time semilog plot. For a given stress load, B and C are 
unknown and may be obtained as before. 

The discussion thus far has centered on a life test where the sample was 
instantaneously loaded to the target tensile stress at some time to. In reality this is not 
possible as there must be a loading rate associated with applying stress to the desired level. 
Thus our simplest physical load history becomes a constant loading rate as depicted in 
Figure A-2. For a loading rate L beginning at time to and ending at tf the stress and 
intrinsic time for the power law breakdown rule are shown to be: 





(A.ll) 



Producing a reliability of, 




(A. 12) 




a ■ 



yields a t «a(b+l) and 



_ ?(p + l) 

P/ “ * P 




(A.13) 
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t 

FIGURE A-2: RAMP LOADING 

These parameters are still in terms of time or life. To put this relationship in terms of 
strength we note that o ult - Lt f for a sample realized in strength, which when substituted 

and compared with a reliability in terms of stress, 




R(o) - exp 



-(f)' 



which for the power law produces a location parameter of 

p, -[i?(p + l)A p p 

Similarly, the exponential law breakdown rule shows that 




or in terms of stress 



R(o ) - exp 




(A.14) 



(A.15) 



(A.16) 



(A.17) 



(A- 18) 
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which produces a shape parameter a=a and location parameter 



P, 



LtB o 




(A. 19) 



Similar derivations may be performed for a variety of break down rules. Difficulties 
may arise in interpreting the break down rule parameters as the convenience of linearity 
may not always be present. Also the application of various load histories was kept to the 
basic stress rupture and constant load rate, but any such history may be applied. 
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APPENDIX B: SERVICE LIFE PREDICTION 



The utility of the discussed relations between strength, life, fiber, and composite may 
be realized upon application to an operational structure. The proposed increase of the P-3 
ZFW was discussed in Reference 1, which provided a finite element analysis of the 
stresses involved. The wing box was marked as a high stress area which would need 
reinforcement should the ZFW be increased. The feasibility of replacing the P-3 wing box, 
shown in Figure 2.1, with composite materials must consider the impact on service life of 
the aircraft. To perform such a study, the techniques discussed in this work will be applied 
in a manner which will hopefully bring together the material covered. 

The basis for the analysis technique discussed in Chapter IV was a commonalty of 
load histories for both fiber-fiber and fiber -composite comparisons. The stress rupture 
load history was chosen as a base for it's relative simplicity and ease of interpretation. 
Once again the power law form of the breakdown rule k will used, emphasizing that a 

correct form must be established experimentally. Note that the stress history experienced 
by the composite structure will affect the degraded time x x through a convolution with the 

chosen k, as shown by: 



Given k and o x (t)., Eq. B.l may be solved for the degraded time for a given load 
history. A larger stress level will provide a larger z x , corresponding to an increase in the 
ZFW. To find the equivalent stress rupture load, set x x - x sr and solve the integrand (Eq. 
B.2) for the equivalent time t x at the desired stress rupture load level. This time t x is the 
equivalent life at the chosen stress rupture load o sr (t), for the service load history o x (t). 




(B.l) 
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(B.2) 



X x m t S r “ 

An idealized load history for the P-3 wing box may be represented by Figure B.l. 
The wing box experiences a gradual stress increase during takeoff as the aircraft weight is 
shifted from the landing gear to the wing generated lift. As fuel is depleted the stress 
steadily decreases during the flight, until landing. The stress spike experienced during 
landing may be significant, and is dependent upon the zero fuel weight. Also, variations in 
the flight stress level can be expected from turbulence, operational maneuvering, and other 
changes in flight conditions. An increase in the ZFW would have the effect of shifting the 
entire stress curve to a higher stress level. What affect this will have on the service life 
must be determined through probabilistic modeling. 




FIGURE B.l: TYPICAL FLIGHT LOAD CYCLE OF P-3 WING BOX 



The cumulation of these flight load cycles represents the total degraded life attributable 
for a load history. To calculate the equivalent stress rupture time for a single flight, the 
flight loading of Figure B.l may be modeled as a series of constant load rates. The 
reduced time per flight is determined by solving: 
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Which may be solved using the power law form of K derived for a constant load rate as Eq. 
A.10. Eq. B.3 reduces to: 



1 


iu\ p (h->or' 


^ p ('3-'2r' + | 




t 


1-4/ P + 1 1-4/ p + 1 


\) p+i 1 


(a) p + l 



(B.4) 

x x in turn is equated to the stress rupture form of Eq. A.3, which provides an equivalent 
stress rupture failure time for an arbitrary stress level. 







(B.5) 



Thus an equivalent stress rupture life may be determined for the specified service load 
history. Fluctuations in the load history may be substituted into Eq. B.l and examined 
similarly. If an appropriate model is known for a composite structure, the problem may be 
worked in reverse to predict the effect of changes in stress on service life. For the case of 
the composite P-3 wing box, knowing the fiber statistics and strength-life model would 
allow for prediction of the effect of increasing the ZFW on service life. 
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APPENDIX C: LIKELIHOOD AND CALCULATION OF 
MAXIMUM LIKELIHOOD ESTIMATORS 



A. BACKGROUND 

Likelihood is defined as the joint probability density function for the random variables 
xi with parameters 0, 

L(X;0) « FT f(x,;6) 

if (C.l) 

where f(x/;0) is the probability density function (pdf) and n the number of samples in the 
data set [12]. The value 0 for which f(jc/;0) is a maximum is called the maximum 
likelihood estimate (MLE) of 0. When applied to the two parameter Weibull distribution 
with shape parameter a and location parameter fl, the likelihood function takes on the form, 



L(*;a > P)-fj/(*,;a,|3) 



where the pdf is defined by 



(C.2) 



/(*,;<*. P) 




(C.3) 



The likelihood defines the probability that the chosen parameters are in fact the correct 
parameters which describe the distribution. The MLE of the pdf indicates the most "likely" 
value of the parameters for the chosen sample. If the partial derivatives of L(Af;a,P) exist, 
then the MLE is equivalent to the solution of the equations 



0L 

da 



3L 

0(3 



0 



For the Weibull distribution the MLE is calculated to be 



(C.4) 
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The solution of Eq. (C.4), depending on the chosen distribution and the number of 
parameters, can be quite tedious and may not produce a closed form solution. If censoring 
of the data is performed, the likelihood function becomes even more cumbersome. As an 
alternative to performing the calculation of Eq. C.5 a graphical or numerical approach may 
be used. This method produces the desired MLE in a fashion more easily interpreted while 
giving insight as to the effects of data censoring. Note that visual interpretation is limited to 
one or two parameter distributions. 

To produce a visual representation of the MLE, values of the likelihood function are 
calculated for a desired range of parameters. This method requires a general knowledge of 
the magnitude of the parameters involved. Once the likelihood is calculated, a contour and 
three dimensional representation may be produced. An example of this procedure is given 
as Figures C.l and C.2. If a great deal of censoring is performed, the likelihood will not 
be as easily interpreted. This is where a visualization of the likelihood surface is helpful in 
determining the most likely values for the parameters. The effects of censoring are 
discussed in Chapter II. 

To investigate the effects of varying parameters and censoring, a set of programs 
written in MATLAB is provided. The program MLE takes a set of exact or censored data 
and determines the MLE parameters for a Weibull distribution. The accuracy of the 
estimated parameters is dependent upon the chosen increment for the shape and location 
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parameters. To simulate a set of Weibull data use the program SIMDATA, with inputs of 
sample size, and shape and location parameters. A flowchart of the software setup is given 
in Figure C.3. 



Realized Data o, Underlying Parameter - 



Marginal Likelihood for F(beta) 




Simulated Data Set 





Likelihood Contour Levels 




alpha (normalized) 

FIGURE C.l MLE CALCULATION OF A WEIBULL DISTRIBUTION, 
N=64, a=5, p=20. 
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Likelihood Contour Levels 
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Likelihood Surface 
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B. SOFTWARE 
% PROGRAM MLE 

% This program performs a Maximum Likelihood Estimator analysis for a given set of data 
% points. This version should be used for analysis of data with small shape parameters 
% which cause a large data spread over a logarithmic scale. A Weibull distribution is 
% assumed and the realized data is ordered using the Expected Rank method. Data is read 
% from data files for exact, right, and left censored data named xi, xr, and xl respectively. 
% If known, population data may be stored as X and Fstar. Program 'Simdata' may be 
% run prior to this program to produce simulated data sets. Ensure the appropriate data file 
% contains a zero if there is no such censoring (i.e. xr=0 if no right censor) 

load xi; % Exact (realized) data points. 

load xr; % Right-censored data points. 

load xv; % Right-censored, Type I-v data points. 

load xl; % Left-censored data points. 

load X; % Population for simulated data sets. 

load Fstarx; % Underlying F* values for simulated data sets. 

load Fstarc; % Censored data F* values 

C=input('Enter type of censor(none=0,interval=l,right=2,left=3): '); 

V=input('If Type I-v present Enter 1: '); 



if V==l 
xr=[xr,xv]; 
end 



if C==0 

m=length(xi); % Number of realized data points. 
n=m; 

elseif C==l 

r=length(xr); % Number of right-censored data points. 
l=length(xl); % Number of left-censored data points. 
m=length(xi); % Number of realized data points. 
n=m+l+r; 
elseif C==2 
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m=length(xi); % Number of realized data points. 
r=length(xr); % Number of right-censored data points. 
n=m+r; 
elseif C==3 

m=length(xi); % Number of realized data points. 
l=length(xl); % Number of left-censored data points. 
n=m+l; 
end 

% Compute the expected rank for each of the samples, and transform to 
% the ln-ln space F* to linearize CDF. Note that the expected rank 
% depends only on the number of samples in the data set. 

erank=[l:n]/(n+l); % Expected rank of sample. 

Fstarex=(log(-log(l-erank)))'; % Expected sample F* values, 
if C==0 
Fstar=Fstarex; 
elseif C==l 

Fstar=Fstarex(l+l:n-r); % Interval data set F* values, 
elseif C==2 

Fstar=Fstarex(l:m); % Right-censored data set F* values, 

elseif C==3 

Fstar=Fstarex(l+l:n); % Left-censored data set F* values, 
end 

% The likelihood calculations performed below produce a square matrix of size 
% 'step.' Rows and columns are determined using increments of alpha and 
% beta respectively. Thus, to obtain a good estimate of the MLE one must have 
% an approximate range for alpha and beta. If not, the program will have to be 
% executed several times before the correct peak can be focused on. 

echo on 

% Input the lower and upper bounds of alpha and beta and desired increment, 
echo off 
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% al=lower alpha, a2=upper alpha, da=alpha increment 
% bl=lower beta, b2=upper beta, db=beta increment 
% a=alpha 
% b=beta 

al=input('Input the lower bound for alpha: '); 
a2=input('Input the upper bound for alpha: '); 
bl=input('Input the lower bound for beta: '); 
b2=input('Input the upper bound for beta: '); 
step=input('Input desired number of steps for alpha and beta: '); 
da=(a2-al)/step; 
db=(b2-bl)/step; 

% Determine parameter increment used and allocate memory space for matrices. 

na=l+(a2-al)/da; 

nb=l+(b2-bl)/db; 

L=zeros(na,nb); 

FA=zeros(l,na); 

FB=zeros(l,nb); 

alpha=al:da:a2; 

beta=bl:db:b2; 

% Compute the likelihood function 



if C==0 

Lex=lnLexact(a l,a2,da,b 1 ,b2,db,xi,m); 
L=Lex; 
elseif C==l 

Lex=lnLexact(a 1 ,a2,da,bl ,b2,db,xi,m); 
Lr=lnLright(al,a2,da,bl,b2,db,xr,r); 
Ll=lnLleft(al,a2,da,bl,b2,db,xl,l); 
L=Lex+Lr+Ll; 
elseif C==2 
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Lex=lnLexact(al,a2,da,bl,b2,db,xi,m); 

Lr=lnLright(al,a2,da,bl,b2,db,xr,r); 

L=Lex+Lr; 
elseif C==3 

Lex=lnLexact(a 1 ,a2,da,b l,b2,db,xi,m); 

Ll=lnLleft(al,a2,da,bl,b2,db,xl,l); 

L=Lex+Ll; 

end 

% Determine the volume under the likelihood surface 
vol=0; 
for i=l:step 
for j=l:step 

avgL=(L(i,j)+L(i+l > j)+L(i,j+l)+L(i+l,j+l))/4; 
vol=vol +(da * db * avgL); 
end 
end 

%L=exp(L)./(-vol); 

L=exp(L); 

% Rearrange likelihood array so (1,1) is plotted at the origin on a contour plot 
% (MATLAB plots (1,1) element in the upper left comer of a contour, not at 
% the origin as one would expect, not needed with version 4.0) 
i=0; 

for k=na:-l:l 
j=0; 
i=i+l; 

for b=bl:db:b2 

j=j+i; 

Lp(k,j)=L(i,j); 

end 

end 
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% Compute the marginal likelihood function, F(alpha) and the approximate 
% area under the curve. 

areafa=0.0; 
for i=l:na 

fa(i)=sum(L(i, :)); 
areafa=areafa+fa(i)*da; 
end 

fa=fa/areafa; 

[fax,i]=max(fa); 

alphamax=alpha(i) 

% Compute the marginal likelihood function, F(beta) 

% and the approximate area under the curve 

dbl=(logl 0(b2)-log 1 0(b 1 ))/step; 

areafb=0.0; 
for j=l:nb 

fb(j)=sum(L(:,j)); 

areafb=areafb+fb(j)*dbl; 

end 

fb=fb/areafb; 

[fbx,i]=max(fb); 

betamax=beta(i) 

betalog=loglO(beta); 

% If simulated data sets have population values of X and F* stored in 
% data files X and Fstarx respectively, then Enter 1 below to add 
% a plot of the true rank. Program Simdata will create these data 
% files automatically for simulations. 

Q=input('Enter 1 if true ranking data is available, else hit Enter key: '); 
clg 
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if 0=1 
subplot(221) 

semilogx(X,Fstarx,'-’,xi,Fstar,'o'),grid 

title('Realized Data (o) and Underlying Parameters (-) in Linearized Domain') 
xlabel('Simulated Data Set'), ylabel('F* Domain') 
else 

subplot(221) 

semilogx(xi,Fstar,'o'),grid 
titleCRealized Data in Linearized Domain') 
xlabel('Data Set'), ylabel('F* Domain') 
end 

subplot(222) 

plot(betalog,fb),grid,title('Normalized Marginal Likelihood for F(beta)') 
xlabel('log beta'), ylabel('F(beta)') 

subplot(223) 

plot(alpha,fa),grid,title('Normalized Marginal Likelihood for F(alpha)') 
xlabel('alpha'), ylabel('F(alpha)') 

subplot(224) 

contour(Lp, 10, beta, alpha), grid, title('Likelihood Contour Levels') 

xlabelCbeta'), ylabel('alpha') 

pause 

dg 

subplot(211) 

contour(Lp, 10, beta, alpha), grid, title('Likelihood Contour Levels') 
xlabelCbeta'), ylabel('alpha') 

subplot(212) 

mesh(Lp) 

title('Likelihood Surface for Selected alpha-beta Range') 
pause 
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function[lnL]=lnLexact(al,a2,da,bl,b2,db,xi,m) 



% Compute the likelihood function for realized data points 
i=0; 

for a=al:da:a2 
j=0; 
i=i+l; 

forb=bl:db:b2; 

j=j+i; 

lnL(i,j)=m*log(a/b)+sum((a-l).*log(xi./b)-(xi./b). A a); 

end 

end 

function[lnLl]=lnLleft(al,a2,da,bl,b2,db,xl,l) 

% Compute the likelihood function for left censored data. Inputs include lower and upper 
% alpha and beta and their respective increments, xl is the value of data points which have 
% been censored to the left with 1 equal to the number of censored points. 

i=0; 

for a=al:da:a2 
j=0; 
i=i+l; 

for b=bl:db:b2; 

j=j+i; 

lnF=log(l-exp(-(xl./b). A a)); 

lnLl(i,j)=sum(lnF); 

end 

end 

function[lnLr]=Lright(al,a2,da,bl,b2,db,xr,r) 

% Compute the likelihood function for right censored data. Inputs include lower and upper 
% alpha and beta and their respective increments, xr is the value of all non-realized data 
% points which have been censored with r equal to the number of censored data points. 



63 



i=0; 

for a=al:da:a2 
j=0; 
i=i+l; 

for b=bl:db:b2; 

j=j+i; 

lnLr(i,j)=sum(-(xr./b). A a); 

end 

end 



% PROGRAM SIMDATA 

% This program produces simulated data which can be used for analysis of of a Weibull 
% distribution model. Inputs include the desired population size and the underlying shape 
% and location parameters, N, alpha, and beta respectively. Outputs include the following 
% column vectors: 

% X = [population] 

% xi = [exact data] 

% xt = [right censored data] 

% xv = [Type I-v censored data] 

% xl = [left censored data] 

% Fx = [True rank] 

% Fstarx = [Underlying F* ranking for simulated data] 

% Fstar = [Expected ranking F* values] 

clear 

clg,hold off 

N=input('Enter the population size: '); 
a=input('Enter desired underlying alpha: '); 
b=input('Enter desired underlying beta: '); 

% Simulate data population and plot F* for expected and true rank vs population 
[X,Fx,Fstar,Fstarx]=population(N,a,b); 
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% Function to produce censored data from expected ranking 
[xi,xr,xv,xl,Fstarc]=censor(X,N,Fstar, Fstarx); 
hold off 

keep=input('Enter 1 to save data: '); 
if keep==l 
save X X /ascii 
save xi xi /ascii 
save xr xr /ascii 
save xv xv /ascii 
save xl xl /ascii 
save Fx Fx /ascii 
save Fstarx Fstarx /ascii 
save Fstar Fstar /ascii 
save Fstarc Fstarc /ascii 
end 

function[X,Fx,Fstar,Fstarx]=population(n,alpha,beta) 

% This function produces a population N characterized by a Weibull distn with shape and 
% location parameters alpha and beta. A plot of the true and expected rank is produced for 
% comparison purposes. 

% Generate a set of random numbers and assign them as the probability of 
% failure F(xi). Order the probabilities F(xi) and transform into the 
% F* space. 

Fx=rand(l,n); 

Fx=(sort(Fx))'; 

Fstarx=log(-log( 1 -Fx)); 

% Solve for the values of X corresponding to the underlying rank 
% given desired shape and location parameters. 

X=beta*(-log(l-Fx)). A (l/alpha); 
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% Compute the expected rank for each of the samples, and transform to the ln-ln space F* 
%to linearize the CDF. Note that the expected rank depends only on the number of 
% samples in the data set. 

erank=[l:n]/(n+l); 

Fstar=log(-log(l -erank))'; 

% Plot the results of the expected rank F and the true rank values F(xi) 

% versus the log(xi) in the linearized F* space. 

semilogx(X,Fstar,'o',X,Fstarx,'-'),grid 
xlabel('Simulated Data Set'), ylabel('F* Domain') 
pause 
end 

functioned ,xr,xv,xl,Fstarc]=censor(X,N,Fstar,Fstarx) 

% This function allows input of a censoring scheme given a simulated 
% population. Right, left, and interval censoring are considered. 

% Type I (or time) censoring occurs if a censor time is prespecified. 

% Type I-v occurs if some outside influence causes unexpected failure. 

% Type II (or failure) censoring occurs if the test is stopped after the 
% rth failure. 

right=input('Enter 1 for right censor data input: '); 
if right==l 

q=input('Enter 1 for Type I censoring or 2 for Type II censoring: '); 
if q==l % Type I censoring 

t=input('Enter desired value of highest realized data point: '); 
for k=l:N 
ifX(k)<=t 

xi(k)=X(k); % xi=realized random variable 
end 
end 
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xi=xi'; 

m=length(xi) ; % Number of realized data points 

if m==N % No data censored 

Fstarc=Fstar; xi=X; 
elseif m==0 % All data censored 

xr=X; censor=l; 
else 

r=N-m; % Number of right censored data points 

Fstarr=Fstar(l:m)'; % Non-censored F* 

Fstarc=Fstarr; 
if r==0, xr=zeros(l,l); 

else, xr=ones(r,l); 
end 

xr=xi(m)*xr; % Assign highest RRV to right censored data 
end 

elseif q==2 % Type II censoring 
m=input('Enter number of realized data points: '); 
if m==N % No data censored 

Fstarc=Fstar; xi=X; 
elseif m==0, 

xr=X; censor=l; % All data censored 
else 

r=N-m; % Number of right censored data points 

Fstarr=Fstar(l:m)'; % Non-censored F* 

Fstarc=Fstarr; 
if r==0, xr=zeros(l,l); 

else, xr=zeros(r,l); 
end 

xi=X(l:m); 
for i=l:r 

xr(i)=xi(m); % Assign highest RRV to right censored data 
end 
end 
end 
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else right=0; 
end 



if xr==0 
Fstarc=Fstar; 
xi=X; 

elseif right==l 
Fstarc=Fstarr; 
end 

left=input('Enter 1 for left censor data input: '); 
if left==l 

l=input('Enter number of left censored data points: '); 
if 1==N, % All data censored on left 

xl=xi; censor=l; 

elseif 1==0 % No data censored 

Fstarc=Fstar; 
else 

Fstarl=Fstar(l+l:N)'; 

xl=ones(l,l); 

xi=X(l+l:N); 

xl=xi(l)*xl; 

Fstarc=Fstarl; 

end 

else left=0; 
end 

% Determine if interval is present 
if left==l 
if right ==1 
intervals 1; 

Fstarint=Fstar(l+ 1 :m); 
xi=X(l+l:m); 

Fstarc=Fstarint; 
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end 

end 



% Determine if no censoring is desired 
if left==0 
if right==0 
Fstarc=Fstar; 
xi=X; 
end 
end 

v=input('Enter 1 if Type I-v censoring desired: '); 
if v==l 

v=input('Enter number of censor points: '); 
for i=l:v 

xv(i)=input('Enter value at censor point: '); 
end 

c=length(Fstarc); 

x=length(xi); 

Fstarc=Fstarc(l:c-v); 

xi=xi(l:x-v); 

censor=0; 

end 

clg, hold off 
if censor==l 

semilogx(X,Fstar,'o',X,Fstarx,'-'),grid 

else 

Fstarc; 

semilogx(X,Fstar,'o',X,Fstarx,'-',xi,Fstarc,V),grid 

end 

title('F* for Expected Rank (o), True Rank (-), and Realized (+) Data') 
xlabel('Simulated Data Set'), ylabel('F* Domain') 
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APPENDIX D: FIBER LOADING PROCEDURES 



Note: The following is an updated version of the loading procedures found in 
Reference 6. A complete mock up of the experimental facility is available for practicing this 
procedure. Its use is highly recommended and significantly reduces errors which produce 
accidental Type I-v censoring of fiber data. 

A. INITIAL SETUP 

1. Ensure power supplies for load cell, elevator, HP-85, and HP-3497A Data 
Acquisition Unit are turned on at least 24 hours prior to loading fibers to help counter 
hysteresis effects. The load cell power supply should be set at approximately 7.5 volts 
using the HP-3497A as a digital readout. 

2. Ensure there are no bubbles present in the hydraulic line connecting the two 
syringes. If so disconnect, empty, and refill with water. This ensures a more constant 
loading rate among fibers. Check that the loader operates smoothly and is capable of 

3 . Check integrity of load cell. Load platform should be very carefully screwed into 
the top of the cell without over tightening. These load cells are easily damaged so be 
careful. 

4. Verify weight of vials to be loaded using balance to nearest 1/10000 gm. 

5. Turn on HP plotter and load paper. Ensure pens are loaded and in working 
condition. Also check paper supply for HP-85 printer. 

6. Adjust the hot filament wire power supply by burning through several discarded 
mounts. The wire should glow a dull red and easily bum through the paper mount 
without causing a fire. Turn the power supply off until needed. 

7 . Prepare work area by removing all unnecessary materials. 
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B LOADER STARTUP AND CALIBRATION 



Overview: The calibration procedure uses a program found in Reference 10 which 
determines the slop and intercept of the calibration curve for the 150 gram load cell. 
Calibration weights are used in even steps to obtain a linear plot of the calibration curve. 
Calibrated weights ranging from 5 to 25 grams are normally used in 5 gram increments. 

1 . Insert tape labeled BELL FIBER TEST into HP-85. 

2. Type CLEAR to clear the screen. 

3. Type LOAD "LDCALB" 

4. Type RUN or press RUN key. 

5 . Answer the prompts on the screen. 

6. Enter the number of calibration points required (usually 5). 

7 . Enter load level to be tested. 

8 . Place calibration weight on center of load cell. 

A. Press ENTER after weight placement. 

B. When system stop reading data remove weight. 

C. Repeat steps 7 and 8 for each calibration point. 

9 . Enter plot axes data. 

A. Enter maximum load to be plotted on calibration curve. 

B. Enter minimum load to be plotted on calibration curve. 

C. Enter maximum X value (read from computer printout, approximately -1.0). 

D. Enter minimum X value (read from computer printout, approximately -4.0). 

10. If full size plot is desired, set up plotter with paper and follow cues on screen. 

1 1 . Repeat the load calibration procedure several times and compare the calibration 
coefficients obtained looking for consistency. The slope coefficient will be used in the 
following procedure. 



C. FILAMENT LOADING 

Overview: The program LOAD5 is designed to measure and record an incremental 
load applied to an unloaded fiber. The desired tensile load is slowly added to the fiber and 
recorded over the 30 second loading period. Fibers which fail during this loading period 
will have the failure load recorded. Fibers which survive the loading process will be 
registered by the Fiber Monitoring System and continue with the life test. 

1 . Insert tape labeled LOADER MOD 1 into HP-85. 

A. Type CLEAR and press ENDLINE key. 

B. Enter LOAD "LOAD5". 

C. Type 251 W=(desired load plus 0.5 gm). This will adjust the y-axis during the 
plotting phase. 

D. Type 55 A=(value of slope coefficient obtained from calibration). 

E. Type RUN 

2 . Carefully remove Plexiglas cover from sample rack. 

3 . Without disturbing adjacent fibers (above, below, and to the side), lower stabilizer 
bolts as far as possible on all stations to be loaded. Loosen the top set screw only 
leaving the bottom screw as a guide for returning the stabilizer bolt to it's original 
position. 

4 . Remove the dead weight or paper block from the station being loaded. 

5 . Adjust the spring so the flag will operate the optical trigger smoothly. Note that the 
Data Acquisition System registers movement of the flag. 

6 . Using tweezers, grasp the fiber to be loaded by one end of it's mount and hang on 
the hook hanging from the optical flag. 

7 . Carefully hang the weight vial from the mount bottom. 

8. Verify the optical flag does not rub against the trigger. 
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9. Adjust the load cell so that the platform is clear from both the hanging vial and the 
stabilizer bolt. 

10. Select function key K1 (Date) on the HP-85 and enter data of the form ddmmmyy. 

11. Select function key K2 (ADJ B) on the HP-85 to adjust the load cell bias. Note the 
weight registered by the computer. If the load is < 0.003 grams, proceed to the next 
step, else perform the bias adjustment until load is within tolerances. 

12. Raise the loading platform such that the vial is situated on the center of the platform 
and the compliance spring and fiber mount cany no load. 

13. Verify that the Monitoring System record the station as being unloaded. 

14. Ensure the direction switch on the loader control panel is in the DOWN POSITION 
and turn the control power switch OFF. 

15. Select function key K3 (LOAD) on the HP-85 and follow cues on screen. 

A. Type in the sample number. 

B. Type in the station number. 

C. Press CONT and verify the weight of the vial shown by the HP-85 is correct for 
the station being loaded. 

16. Turn the power supply for the filament wire ON. 

17. Without touching or even approaching the fiber to be loaded, bum through both 
sides of the fiber mount at the bottom of the slot. A thin wooden dowel or similar device 
may help stabilize the mount from spinning. 

18. Switch the hot wire filament power supply OFF and return to it's holder. 

19. Simultaneously press CONT and switch the loader control power ON to continue 
the program and start a 2 minute timer. 

20. If fiber breaks during loading (If fiber survives loading proceed to step 21): 

A. Carefully pull fiber weight from load cell. 

B. Turn loader power switch to OFF. 
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C. Replace weight on load cell. 

D. Continue running program to 1 minute 45 seconds. 

E. Select PAUSE on the HP-85. 

F. Type CONT 100 and press END LINE on HP-85. 

G. Select function key K4 (PLOT) while ensuring that the error light on the plotter is 
off. If the error light is on, recycle the power switch. 

H. Follow cues on screen to operate the plotter. 

I. Select PAUSE. 

J. Ensure the Fiber Monitoring System records the appropriate results and record 
fiber load data in the Fiber Life Data Book. 

K. Type RUN and continue from step 4 for additional fiber loading, else proceed to 
step 22. 

21 . If the fiber survives the loading process: 

A. Verify the Monitoring System records the station as being loaded when the flag 
blocks the optical trigger. 

B. When the fiber is completely loaded the HP-85 will beep twice. At this time turn 
the loader control switch OFF 

C. Select PAUSE on the HP-85. 

D. Type CONT 100 on the HP-85. 

E. Select function key K4 (PLOT) to plot the loading process while ensuring that the 
error light on the plotter is off. If the error light is on, recycle the power switch. 

F. Follow cues on screen to operate plotter. 

G. Select PAUSE. 

H. Ensure the Fiber Monitoring System records the appropriate results and record 
fiber load data in the Fiber Life Data Book. 
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I. Type RUN and continue from step 4 for additional fiber loading, else proceed to 
step 22. 

22. Adjust stabilizer bolts such that there is minimal spacing between the vial weight 
and the top of the bolt. 

23. Replace Plexiglas covers and power down equipment. 

24. Verify the Fiber Monitoring System is reporting the correct status for the stations 
loaded. 
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APPENDIX E: LIFE MONITOR SYSTEM DESCRIPTION 



Fiber life data was collected using a data acquisition system which monitors 512 
suspended fibers at once. Driven by a PC-AT computer, the acquisition software marks 
failure times by monitoring an infrared switch at each station. The switch is closed by a 
flag attached to a spring at one end and the test fiber at the other. The flag cycles up when 
the fiber fails to close the switch. The current time and loading time are then compared to 
produce a life span for the fiber. All fiber data is then manually recorded in spreadsheet 
form. A more detailed description of system components may be found in Reference 10. 

Experience has shown that the system is easy to maintain provided that all notation in 
the Load Manual and Life Data spreadsheet are thoroughly understood. New personnel are 
encouraged to complete a thorough check and update of the acquisition system and related 
data base. Periodic power-up and power-down of the system may be required for various 
reasons. To facilitate this process, the following procedure is recommended. 

A. POWER DOWN ( AC POWER AND BATTERY MODE) 

1. Monitor - OFF 

2. Printer - OFF 

3. Computer - OFF 

4. Optical Switch Interface - OFF 

5. Power Supply - OFF 

B . POWER UP ( AC POWER AVAILABLE ) 

1. Power Supply - ON 

2. Optical Switch Interface - ON 

3. Printer - ON 
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4. Monitor - ON 



5. Computer - ON 

6. Reboot system by pressing CTL, ALT, DELETE together 

C . POWER UP ( BATTERY POWER) 

1. Power Supply - ON 

2. Optical Switch Interface - ON 

3. Printer - ON 

4. Computer - ON 

5. Verify start time on print-out 

6. Reboot system by pressing CTL, ALT, DELETE together 

D . DATA RECORDING PROCEDURE 

1. Fiber failure written to printer, screen, floppy. 

2. Verify physical location of station and confirm failure. 

3. Remove printout and compare with screen data. 

4. Compare fiber data and load information with Load Manual. 

5. Deactivate station with broken sample. 

6. Remove floppy and back up on different computer. 

7. Replace floppy and reboot computer. 

8. Record failure on spreadsheet titled "Life Data mm/dd/yy.spd." 

9. Verify spreadsheet calculated life with monitor system time. 

10. Save updated spreadsheet as a new file using current date. 

11. Delete the oldest file ( a total of five versions should be present). 
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APPENDIX F: TRANSFORMATION SOFTWARE 



The following software package was designed to analyze a chosen form of the break 
down rule. It is currently configured for the power law breakdown rule, however any 
form may be used, so long as all parameters are accounted for. The program was written 
in MATLAB and structured as shown in Figure F.l. 
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Ramp-to-strength-to-life 
Mapping Program 



Load xi (RRV), 
ramp loading data, 
file 



C 



Powerlaw Analysis 
Subroutine 



i 



<|e 



Input Load Rate, slope 
of strength-life 
degradation, and 
intrinsic time 



Determine sample size 

Input estimated shape 
and location 
parameters for data set j 



i 



Determine intrisic strength of each 
RRV based on load rate and 
reference time (deterministic 
approach) 



i 



Output to breakdown rule subroutine 


i 


i 


/ Plot strength-life / 

/ degradation of mapped / 
/ data / 


\ 





i 



Determine the location parameter for 
the stress axis and the shape 
parameter in time. 



Input stress levels for 
strength-to-life 
transformation 



Output for 3D CDF surface based 
on chosen breakdown rule 



x 

C M ) 
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Map the intrinsic strength data to specified stress 
levels. Use the location parameter on the stress axis 
to determine the location parameter in time for each 
stress level. 



r — Ls 

I Return J 



FIGURE F.1A PROGRAM STRENGTHLIFE FLOWCHART 
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FIGURE F.1B PROGRAM STRENGTHLIFE FLOWCHART 
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% PROGRAM STRENGTHLIFE 

% This program allows for a testing of the power form of the breakdown rule. 
% This program will take data for a ramp load as input, determine the required 
% parameters and transfer the data to the stress axis at the desired intrinsic 
% time and to any desired load level. 

% Input data may take the form of simulated or actual data, but all forms must 

% be listed in column vector format and labeled as indicated below. 

clear 

load xi; % RRV from true ranking 

xramp=xi; 

n=length(xramp); 

% alphar = shape parameter for ramp loading condition 
% Br = location parameter for ramp loading condition 

alphar = input('Enter ramp data estimated shape parameter (alphar): '); 

Ur = input('Enter ramp data estimated location parameter (betar): '); 

% Transfer ramp load data to stress axis and life axis using break down rule 

pwrlaw % Power law form of break down rule used 

% Produce load level vectors for tf plot 

Q=ones(n,l); 

sigl=stress(l).*Q; 

sig2=stress(2).*Q; 

sig3=stress(3).*Q; 

sig4=stress(4).*Q; 

sig5=stress(5).*Q; 

thatv=that.*Q; 

time 1 = [tf 1 ( 1 ),tf2( 1 ),tf3( l),tf4( 1 ),tf5( l),that] ; 
time2=[tfl(n),tf2(n),tf3(n),tf4(n),tf5(n),that]; 
sigmal=[stress,A(l)]; 
sigma2=[stress,A(n)]; 
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% Plot stress-life degradation on loglog plot 
clg, hold off 

loglog(tfl,sigl,'o',tf2,sig2,'o’,tf3,sig3,'o',tf4,sig4,'o , ,tf5,sig5,'o',... 
thatv. A, 'x',timel,sigmal,'-',time2,sigma2,'-') 

%title('Realized Data Plot, "o" realized in time, "x" realized in strength') 
pause 

% Plot CDF surface using 3D plotting routine for specified breakdown 
% rule. SRpwrTD uses power law relationship. 

td=input('If a 3-dimensional CDF plot is desired enter 1: '); 
if td==l 
SRpwrTD 
end 

% Save plot data for easy import to other graphing programs 
keep=input('Enter 1 to save 3D plot data and failure times: '); 
if keep==l 
save cdf Ft /ascii 
save tfl tfl /ascii 
save tf2 tf2 /ascii 
save tf3 tf3 /ascii 
save tf4 tf4 /ascii 
save tf5 tf5 /ascii 
end 



% POWER LAW SUBROUTINE 

% Transfer ramp load data to desired stress axis to obtain a stress 
% rupture loading condition. Then transfer data to any desired stress level. 
% Note that the A vector is equivalent to the intrinsic strength of the sample 
% at time t-hat. The parameter a is interpreted as alphat in the SR domain. 
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% Ldot = loading rate (gm/sec) 

% rho = initial slope of stress-life correlation curve 
% that = intrinsic life denoted as "t hat" 

% A = parameter for power law equation 
% a = hazard function parameter 

% 6s = shape parameter in stress at "t hat", required for 3D plot 

Ldot =input('Enter ramp load rate (gm/sec): '); 
rho =input('Enter initial slope value: '); 
that =input('Enter intrinsic life: '); 

A=xramp. A ((rho+l)/rho)./(that*Ldot*(rho+l)). A (l/rho); 

6s=6r A ((rho+ 1 )/rho)/(that * Ldot * (rho+ 1 )) A ( 1/rho) 

a=alphar/(rho+l); 

alphat=a; 

% Transfer stress axis data to 5 defined load levels using power law and 
% determine failure times. Stress rupture form of the power law is used. 



for i=l:5 

stress(i)=input('Enter stress level for life analysis: '); 
end 

tfl=that.*(A/stress(l)). A rho; 
tf2=that.*(A/stress(2)). A rho; 
tf3=that.*(A/stress(3)). A rho; 
tf4=that.*(A/stress(4)). A rho; 
tf5=that. *(A/stress(5)). A rho; 

6t 1 = that * (6s/stress(l )). A rho 
6t2=that * (6s/stress(2)). A rho 
Bt3=that*(Bs/stress(3)). A rho 
Bt4=that*(Bs/stress(4)). A rho 
Bt5=that * (Bs/stress(5)). A rho 



83 



% STRESS RUPTURE THREE DIMENSIONAL PLOT ROUTINE 
% 3D plot routine for stress rupture loading of a fiber population assuming 
% a Weibull distn. This version uses the power law form of the break 
% down rule assuming a linear relationship between stress and time in a loglog 
% plot. 

% Input Variables 

k = input('Enter desired mesh increment: '); 

% Determine max and min values for the stress and time axis. Note the 
% use of logrithmic increments. Program may be modified to allow for 
% input of stress and time min/max values by swapping comment symbols (%) 

%smax = input('Enter max stress for 3D plot: '); 

%smin = input('Enter min stress for 3D plot: ’); 
smaxl = [max(stress) max(A)]; 
smax - max(smaxl); 
sminl = [min(A) min(stress)]; 
smin = min(sminl); 
ds = (Iogl0(smax)-logl0(smin))/k; 
sexp = [loglO(smin):ds:loglO(smax)]'; 

%that = input('Enter min time for 3D plot: '); 

%tmax = input('Enter max time for 3D plot: '); 
tx = [max(tfl) max(tf2) max(tf3) max(tf4) max(tf5) that]; 
tmax = max(tx)+0.2*max(tx); 
tn = [min(tfl) min(tf2) min(tf3) min(tf4) min(tf5) that]; 
tmin = min(tn); 

dt = (Iogl0(tmax)-logl0(tmin))/k; 

% Determine CDF and PDF for designated stress-time region. For each stress 
% level determine the failure times for each sample. Bm is the location 
% parameter for each stress based on the linear power law form. Note use of 
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% logrithmic increments for si and ti. 



F = zeros(k); 
si = smin; 
for s=l:k 
ti=tmin; 

Bti=that*(Bs/si) A rho; % Location parameter in time for each stress level 
for t=l:k 

Ft(t,s)=l-exp(-(ti/Bti) A alphat); % Weibull form of CDF 
ti=10 A (logl0(ti)+dt); 
end 

si=10 A (logl0(si)+ds); 

end 

clg 

hold off 

M=[-70 15]; %Set perspective for view 
mesh(Ft,M) 

title('Cummulative Distribution Function for Stress Rupture Loading') 
pause 
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